
x = [-1 1 1 -1];
y = [1 1 -1 -1];

vert = [x', y'];
edge = zeros(4, 2);
edge(:, 1) = [1:4]';
edge(:, 2) = mod(edge(:, 1), 4) + 1; 

param = createParam(0.00001, 0.03, 0.0000025, 0.1);

domain2D = createDomain(vert, edge, param);

dim = size(domain2D.vx, 1);

b = assemblyBoundaryVector(domain2D);

res = 50;
[rx, ry] = meshgrid(0:res,0:res);
rx = rx / res * 1.86 - 0.93;
ry = ry / res * 1.86 - 0.93;
rx = reshape(rx, (res+1)*(res+1), 1);
ry = reshape(ry, (res+1)*(res+1), 1);

[A, Z] = assemblyMatrix(domain2D);

theta = A \ b;

u = evaluateInteriorPoint(rx, ry, domain2D, theta);

F = evaluateInteriorDeformationGradient(rx, ry, domain2D, theta);